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Abstract 

We introduce a direct method allowing to solve numerically inverse type problems for 
linear hyperbolic equations. We first consider the reconstruction of the full solution of the 
wave equation posed in x (0, T) - 17 a bounded subset of - from a partial boundary 
observation. We employ a least-squares technique and minimize the L^-norm of the distance 
from the observation to any solution. Taking the hyperbolic equation as the main constraint 
of the problem, the optimality conditions are reduced to a mixed formulation involving both 
the state to reconstruct and a Lagrange multiplier. Under usual geometric optic conditions, 
we show the well-posedness of this mixed formulation (in particular the inf-sup condition) and 
then introduce a numerical approximation based on space-time finite elements discretization. 
We prove the strong convergence of the approximation and then discuss several examples for 
N = 1 and N = 2. The problem of the reconstruction of both the state and the source term 
is also addressed. 


1 Introduction 

Let 17 be a bounded domain of (TV > 1) whose boundary dil is of class and let T > 0. 
We define Qt := 17 x (0,T), Et := x (0,T) and denote by iz = i/{x) the outward unit normal 
to 17 at any point x G 917. We are concerned with inverse type problems for the following linear 
hyperbolic equation 

( ytt-y ■ic{x)Vy)+d{x,t)y = f, {x,t) € Qt 

< y = 0, (a;,i)eET (1) 

I iy{-,0),yt{-,0)) = {yo,yi), a; e 17. 

We assume that c S <17^(17, K) with c{x) > cq > 0 in 17, d S L°°{Qt), {yo, Vi) G V := iLg (17) xL^(17) 
and / G L‘^{Qt)- 

For any {yo,yi) G V and any / G L‘^{Qt), there exists exactly one solution y to ([^, with 
y G C°([0, T]; iJg (17)) n C^([0, T]; L^(17)) (see [HHH]). In the sequel, for simplicity, we shall use 
the following notation: 

Ly ■■= ytt-y ■ icix)Vy) + d{x,t)y. (2) 
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Let now L be any non empty open subset of dO and let Ly := F x (0, T) C Tit- A typical inverse 
problem for Q is the following one : from an observation or measurement of the normal derivative 
yiy,obs in L'^{Tt) on the sub-domain Ft’, we want to recover a solution y of the boundary value 
problem 0 such that its normal derivative coincides with the observation on F^. 

We denote by d^y = Vj/ • n the normal derivative of y on dD. Introducing the operator 
P : Z ^ L'^{Qt) X L'^(Tt) defined by Py := (Ly, c(a;)i9yy|ry) and the Hilbert space Z defined by 
0-0, the problem is reformulated as : 

find y G Z solution of P y = {f,c{x)y^^obs)- (IP) 


From the unique continuation property for 0> if the set F^ satishes some geometric conditions 
and if yv^obs is a restriction to Fy of a normal derivative of a solution of 0, then the problem is 
well-posed in the sense that the state y corresponding to the pair {yu,obs, f) is unique. 

In view of the unavoidable uncertainties on the data yu,obs (coming from measurements, numer¬ 


ical approximations, etc), the problem (IP) needs to be relaxed. In this respect, the most natural 
(and widely used in practice) approach consists in introducing the following extremal problem (of 
least-squares type) 


minimize over V J(yo,yi) := ^llc(x)(di,y - yi,^obs)Hl2(rj,) 
where y solves 0, 


(LSJ 


since y is uniquely and fully determined from / and the data (yo, yi). Here the constraint c{x){d^y— 
yv,obs) = 0 in L‘^{Tt) is relaxed; however, if y^^obs is a restriction to Ft of the normal derivative of 
a solution of 0> then problems ( |F,S[ ) and ([70 coincide. A minimizing sequence for J in V is easily 
defined in term of the solution of an auxiliary adjoint problem. Apart from a possible low decrease 
of the sequence near extrema, the main drawback, when one wants to prove the convergence of a 
discrete approximation is that, it is not in general possible to minimize over a discrete subspace of 
{y; Ly —/ = 0} subject to the equality (in L'^{Qt)) Ly — f = 0. Therefore, a classical trick consists 
first in discretizing the functional J and the system 0; this raises the issue of uniform coercivity 
property (typically here some uniform discrete observability inequality for the adjoint solution) 
of the discrete functional with respect to the approximation parameter. As far as we know, this 
delicate issue has received answers only for specific and somehow academic situations (uniform 
Cartesian approximation of 0, constant coefficients in 0, etc). We refer to [HI [HI EOl [24] and 
the references therein. 


More recently, a different method to solve inverse type problems like (IP) has emerged and use 


the so called Luenberger type observers: roughly, this method consists in defining, from the obser¬ 
vation on Ft, an auxiliary boundary value problem whose solution possesses the same asymptotic 
behavior in time than the solution of 0: the use of the reversibility of the hyperbolic equation 
then allows to reconstruct the initial data (yo,yi)- We refer to [HIES] and the references therein. 
However, for the same reasons, from a numerically point of view, these methods require to prove 
uniform discrete observability properties. 

In a series of works, Klibanov and co-workers use different approaches to solve inverse problems 
(we refer to [19] and the references therein): they advocate in particular the quasi-reversibility 
method which reads as follows : for any e > 0, find y^ G A the solution of 


{Pye,Py)mQTUL^(TT)+ £{ye,y)A = {if,c{x)yu,obs),Py)T2(^Q^)xLUrT) ’ ^ y £ A (QRe) 


where A denotes a Hilbert space subset of L'^{Qt) so that Py G L‘^{Qt) x L^(Ft) for all y G A and 
e > 0 is a Tikhonov like parameter which ensures the well-posedness. We refer for instance to [H] 
where the lateral Cauchy problem for the wave equation with non constant diffusion is addressed 
within this method. Remark that (QRe ) can be viewed as a least-squares problem since the solution 
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Ds minimizes over A the functional y \\Py - (/, c(a;)?/^,obs)|li 2 (Q^)xL 2 (rj,) + ^WuWa- Eventually, 
if 2/i/,obs is the normal derivative of a restriction to Ft of a solution of Q, the corresponding y^ 
converges in L‘^{Qt) toward to the solution of {IP) as e —>■ 0. There, unlike in Problem {LS), the 


unknown is the state variable y itself (as it is natural for elliptic equations) so that any standard 
numerical methods based on a conformal approximation of the space A together with appropriate 
observability inequalities allow to obtain a convergent approximation of the solution. In particular, 
there is no need to prove discrete observability inequalities. We refer to the book [2]. We also 
mention w where a similar technique has been used recently to solve the inverse obstacle problem 
associated to the Laplace equation, which consists in finding an interior obstacle from boundary 
Cauchy data. 

In the spirit of the works [iiiaiini, we explore the direct resolution of the optimality conditions 


associated to the extremal problem [LS), without Tikhonov parameter while keeping y as the 


unknown of the problem. This strategy, which avoids any iterative process, has been successfully 
applied in the close context of the exact controllability of 0 in m and 0111]. The idea is to 
take into account the state constraint Ly — / = 0 with a Lagrange multiplier. This allows to derive 


explicitly the optimality system associated to (LS) in term of an elliptic mixed formulation and 


therefore reformulate the original problem. Well-posedness of such new formulation is related to a 
unique continuation property for the hyperbolic equation 0 - 

From the observation y^^obs, we also address the simultaneous reconstruction problem of the 
source term / and the solution y: 


find {y,f) G Z X L^{Qt) solution of Py = {f,c{x)y^^obs), 


(IPf) 


where (j/,/) solves 0. Without additional assumption on f G L'^{Qt), the pair {y,f) solution 
of \IPfl is not unique: consider for instance a source term / supported in a set which is near 
n X {T}: from the finite propagation of the solution, the source / will not affect the solution y on 
Tt. On the other hand, a result of Yamamoto and Zhang in m asserts that the uniqueness holds 
true if the source takes the form f{x,t) = o'(t)yt{x), where the smooth time part a is given and 
the spatial part /r is a function. 

We adapt in this work the arguments of m where the observation is distributed in Qt. The 
outline is as follow. In Section]^ we consider the least-squares problem 0 and reconstruct the 
solution of the hyperbolic equation from a partial observation localized on a subset Tr of Et- 
For that, in Section 2.1 we associate to © the equivalent mixed formulation ([^ which relies on 
the optimality conditions of the problem. Assuming that Ft satisfies the classical geometric optic 
condition (Hypothesis 1, see we then show the well-posedness of this mixed formulation, in 

particular, we check the Babuska-Brezzi inf-sup condition (see Theorem 2.1). Interestingly, we also 


derive in Section 2.2 an equivalent dual extremal problem, which reduces the determination of the 


state y to the minimization of a strongly elliptic functional with respect to the Lagrange multiplier. 
In Section we use the uniqueness result m Theorem 2.1] of Yamamoto and Zhang and apply 
the same procedure to recover from a partial observation both the state and the spatial part y of 
the source term assumed in Sectionj^is devoted to the numerical approximation, through 

a conformal space-time finite element discretization. The strong convergence of the approximation 
{Vh, yh} is proved as the discretization parameter h tends to zero. In particular, we discuss the 
discrete inf-sup property of the mixed formulation. We present numerical experiments in Section 
[^for fl = (0,1) and B C an check the robustness and convergence of the approximations in 
agreement with the theoretical part. Section [^concludes with some perspectives: in particular, we 
highlight that the parabolic case can be treated in a similar way. 
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2 Recovering the solution from a partial observation: a 
mixed re-formulation of the problem 

In this section, assuming that the initial data (yojyi) G ^ are unknown, we address the inverse 
problem {IP). Without loss of generality, in view of the linearity of the system Q, we assume 
that the source term / is zero. 

We consider the non empty vector space Z defined by 

Z:={y.ye C{[0,T], H^{n)) O C\[0,T], L^O)), Ly G L^Qt)} (3) 

and then first recall that d^y G L‘^{Et) for any y G Z: precisely (see [2T1 Theorem 4.1, Ch 1]), 
there exists a constant Ct > 0 such that the following holds : 

(4) 


||c(a;)9j,?/||^2(r,j,) < Ct (^|| (?/(•, 0), 0))|1^ + j , ^y G Z. 

We then introduce the following hypothesis : 

Hypothesis 1 There exists a constant Cobs = C{uj,T, ||c||(;;i(qj, such that the follow¬ 

ing estimate holds : 


l|y{’) 0), yt(-, 0)11 V < Cobsi l|c(a;)9;,y||^2(r.j,) + ||iy||L2(Q,j,) )> '^V ^ 




Condition is a generalized observability inequality for the solution of the hyperbolic equa¬ 
tion Q. For constant coefficients, this estimate is known to hold if the triplet (r,T, fl) satisfies a 
geometric optic condition. We refer to [1] for the case of constant velocity c. In particular, T must 
be large enough. In the one-dimensional case, for non constant velocity c and potential d, we refer 
to El and the references therein. 

Then, within this hypothesis, for any 77 > 0, we define on Z the bilinear form 


{y,y)z-= c'^{x)doyd^ydadt-\-r] LyLydxdt \/y,y G Z 

J J rjp J J 


(5) 


and we denote the corresponding semi-norm ||y||z := \/(y, y)z- We have the following result : 

Lemma 2.1 Under the hypothesis {%), the space Z is a Hilbert space with the scalar product (•, ■)z 
defined by @- 

PROOF-The two main properties we need to verify are that the semi-norm associated to this inner 
product II • \\z is indeed a norm, and that Z is closed with respect to this norm. The first property 
is a direct consequence of the inequality (H). 

To check the second property, let us consider a convergent sequence {zk}k>i C Z such that 
Zk ^ z in the norm |1 • jj^. We have to see that z G Z. As a consequence of ("H), there exist 
{zq,Zi) G V and / G L'^{Qt) such that ( 2 fc(-, 0), 0)) —>■ (^ 0 , 2 : 1 ) in V and Lzk —>■ / in 

L^(Qt). Therefore, Zk can be considered as a sequence of solutions of the hyperbolic equation 
with convergent initial data and second hand term Lzk —>■ /. 

By the continuous dependence of the solutions of the wave equation on the data, Zk ^ z in 
(7(10, T]; iLg (n)) n T]; L‘^{Ll)), where 2 : is the solution of the hyperbolic equation with initial 

data ( 2 : 0 , 2 : 1 ) G V and second hand term Lz = f G L'^{Qt). Therefore, in view of Q, z G Z. □ 


We consider the following extremal problem : 

J inf J(y) := h\c{x){d^y - y,,,obs)\\l2^rT)^ 


1 


subject to y gW 


(V) 
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where W is the closed subspace of Z defined by 

W:={yGZ; Ly = 0 in L^{Qt)} 
and endowed with the norm of Z. 

The extremal problem © is well posed : the functional J is continuous over W, is strictly 
convex and is such that J(y) —>■ +oo as \\y\\w —> oo. Note also that the solution of © in W does 
not depend on rj. 

We recall that from the definition of Z, Ly belongs to L‘^{Qt)- Furthermore, the uniqueness 
of the solution is lost if the hypothesis is not fulfilled, for instance if T is not large enough. 
Eventually, from (|^, the solution y in Z of (|^ satishes (y(-, 0), yt(-, 0)) G V, so that Problem © 
is equivalent to the minimization of J with respect to iyo^yi) as in problem {IP), Section 1. 


2.1 Direct approach 

In order to solve ©, we have to deal with the constraint equality which appears in the space W. 
Proceeding as in m, we introduce a Lagrange multiplier A G L^{Qt) and the following mixed 
formulation: find (y. A) G Z x LP‘{Qt) solution of 


a{y,y)+ h{y,X) = l{y), 

= 0, 


VyG Z 

yXeL^iQr), 


where 


a : Z X Z 


a{y,y) := // c^{x) d,,y d,,y dadt, 

77 Tt’ 


6 : Z X —>■ // XLydxdt^ 

J J Qt 

/ : ZK, l{y) ■■= // c^{x)y,,^obsd,,y dadt. 
JJtt 


( 6 ) 

(7) 

( 8 ) 

(9) 


System § is nothing else than the optimality system corresponding to the extremal problem ©■ 
Precisely, the following result holds: 

Theorem 2.1 Under the hypothesis 

(i) The mixed formulation © is well-posed. 

(ii) The unique solution (y, A) G Z x L^{Qt) to is the unique saddle-point of the Lagrangian 
C : Z X L‘^{Qt) -G K defined by 

■=\a{y, y) + b{y, A) - l{y). 

(in) We have the estimate 

1 

\\y\\z = (\\c{x)d,,y\\l2^rT)+y\\^y\\h(QT))^ ^ \\cix)yu,obs\\mrT)^ 

I!^IIl2(Qt) <2y'Cn,T + v\\ yu,obs II (Et) ■ 


( 10 ) 


Proof- The proof is based on classical results for saddle point problems (see [1], chapter 4). 
We easily obtain the continuity of the bilinear form a over Z x Z, the continuity of bilinear b 
over Z X L^{Qt) and the continuity of the linear form I over Z. In particular, we get 

II^IU' < \\c(x)y, 2 ,obs\\L^(rT)j l|a||(ZxZ)' < 1, II&II(ZxL2(Qt))' < (H) 
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Moreover, the kernel JC{b) = {y £ Z; b{y, A) = 0 VA S L‘^{Qt)} coincides with W: we easily get 

a{y,y) = \\yfz, Vy e 7V(6) = kh. 

Therefore, in view of [U Theorem 4.2.2], it remains to check the inf-sup constant property : 35 > 0 
such that 

KyA) 


inf sup ■ 


> 5. 


\eL^(QT)yez ||y||z||A||L2(QT) 

We proceed as follows. For any hxed A G L‘^{Qt), we define yo £ Z as the unique solution of 

Lyo = A in Qt, (yo(-,0),yo.t{-,0)) = (0,0) in O, yo = 0 on S-r- 
We get b{yo,X) = ||A||| 2 (q,^) and 

llj/o||| = l|c(a;) 9i,yo|li2(r,j,) + v\\X\\^i^ 2 (^Q^y 
Using Q, the estimate ||c(a;)5^yollL2(rT) < ^C^llA|| l 2 (Qi,) implies that 

b{yo,X) , 1 


( 12 ) 


sup 


> 


> 


yez I| 2 /||z||A||l 2 (q,j,) ||yo||z||A||L2(Qr) + V 


> 0 


leading to the result with S = {Cq„t + v) ■ 

The third point is the consequence of classical estimates (see [4], Theorem 4.2.3.) 

ii2/iu<-ii/iU', mLHQ.)<Ui + —)\\i\\z' 
ao ' 6 \ ao J 

■ r a{y,y) 

ao ini ,, ■ 

yeN{b) llyjj^ 


where 


(13) 


Estimates (111 and the equality ao = 1 to the results. Eventually, from (11), we obtain that 

2 

||A||l2(Qt,) < j\\c{x)y^^obs\\L^{rT) 


and that S > {Cn,T + v) to get (101. □ 

In practice, it is very convenient to ’’augment” the Lagrangian (see m) and consider instead 
the Lagrangian Cr defined for any r > 0 by 

^r{y, A) := ^ar{y, y) V b{y, A) - l{y), 
ar{y,y) := a{y,y) + r\\Ly\\l 2 ^Q^y 

Since 0 ^( 1 /, y) = a{y, y) on W, the Lagrangian £ and £r share the same saddle-point. The positive 
real r is an augmentation parameter. 


Remark 1 Assuming additional hypotheses on the regularity of the solution X, precisely LX G 
L'^{0,T, H~^(Jl)) and (A, At)|t=o,r G L'^{Ll) x we easily check, by writing the optimality 

conditions for £, that the multiplier X satisfies the following relations : 


0 

II 

in 

Qt, 

A = c(x) (5i,y - yu,obs) 

on 

Tt, 

A = 0 

on 

E-r \ Ft, 

0 

II 

-< 

II 

on 

n X {o,r} 


(14) 


Therefore, X (defined in the weak sense) is a exact null controlled solution of the hyperbolic equation 
^ through the boundary control {d„y - obs ) IFt GL^Tt). 
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obs is the normal derivative of a solution of 0 » restricted to Tt, then the unique multi¬ 
plier A must vanish almost everywhere. In that case, we have 


with 


sup inf Cr{y,X) = inf Cr{y,0) = inf Jr{y) 
ag_L2(Qj,) yez yez yez 


1 


2^ W/ y'' iz if ijtz-,uua / \ rj,j ■ 2 ' 

The corresponding variational formulation is then : find y € Z such that 


Jr{y) -z\\c{x){d,,y — yv,obs)\\L'^{TT) 7 )\\^y\\L^(QT)- 


(15) 


ar{y,y)= // c'^{x)d,,yd,,ydadt-\-r Ly Ly dxdt = l{y), My 

JJOt' 


€ z. 


'Ft 

• In the general case, the mixed formulation can he rewritten as follows: find {y, X) G Z x 
L^{Qt) solution of 

j {Pry,Pry)L2(QT)xL^irT) + (-^1/; (Q j,) = {{0, c{x)y,,^obs), PrV) {Qt)x (Ft) > ^ Z, 

\ {Ly,X)L'^{QT) = MXGL‘^{Qt) 

with Pr : Z X LF'(Qt)z L‘^(Qt) defined by P^y := {\/rLy,c{x) d,jy\Y-T)- This formulation may 
be seen as a generalization of the quasi-reversibility formulation \QRe^ , where the variable 
X is adjusted automatically (while the choice of the parameter e in (QRe ) is in general a 
delicate issue). 


System (14) can be used to define a equivalent saddle-point formulation, very suitable at the 


numerical level. Precisely, we define - in view of ( |14[ ) - the space A by 
A := {A : A e C{[0,T];L\n)) n C\[0,T]-, 

LX G = At(-,0) = 0,A|r^ G L^{Tt)}. 


Similarly to Lemma 2.1 we can prove that A is a Hilbert space endowed with the following inner 
product 

(A, A)a := f {LX{t),LX{t))H-^(n)dt-\- ff (?{x)XXdadt, V A, AgA 

t/ 0 J J Y''j’ 

using notably, that the elements of the non-empty vector space A satisfy the inequality 

I1 A||l2(q,j,) < Cn,T\/ (A, A)a (16) 

for some positive constant Cq t which depend on and ||(i||Aoo(Q,j,). In what follows we 

denote ||A||a := ^/{X,X)a- 

Then, for every parameter a G (0,1), we consider the following mixed formulation: 


ar,aiy,y) +baiy,X) = ll,aiy), 
baiy^ X') Cfj(A, A) ^2,aiXf 


vyG z 

V A G A, 


where 


a^^a • Z X Z —^ 
be '. Z X A —y 
Ce : A X A —y 
h,a : T —t 
h,a : A —>■ 


ar',a(2/,y) = (1 - a) // c{x)dvy duy dadt + r // Ly Lydxdt 

J J r 'JP d J 

ba{y,X)= // LyXdxdt — a // c'^{x)d,yyXdadt 

J J (^'j' J J r 'j' 

Ca{X,X)=a f {LX{t),LX{t)) H-^{Q.)dta [[ (?{x)XXdadt 


Ft 


h,a{y) = (1 - a) I c^{x)duyyu,obsdadt 

I Ft 


h.aiX) =-a c{x)y,,nbsXdudt. 


(17) 

(18) 

(19) 

( 20 ) 
( 21 ) 
( 22 ) 
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From the symmetry of Ur^a and Cq, we easily check that this formulation corresponds to the 
saddle point problem : 


{ sup inf/:r.a(y, A), 

Or,a{y, A) := Cr{y, A) — ^||iA||^ 2 (//-i(Qp — -^\\c[x) (A — [d^y — yv,obs)) Ili, 2 (r^)- 


(23) 


Precisely, the following holds true. 


Proposition 2.1 


{17} is well-posed. 


Under the hypothesis for every a S (0,1) the stabilized mixed formulation 
Moreover, the unique pair {y, X) G Z x A satisfies 


■ a 


11 < 


(1 — a)^ + aO I 

9 


\ytJ,obs\\L'^{TT) 


(24) 


with 6 := min (1 — a, r/rf). 


Proof- We easily get the continuity of the bilinear form Or,a, ba and Ca'. 


|ar,a( 2 /,y)| < max - a, ||?/||z|| 2 /||z, Vy, y G Z, 
|&a(2/,a)|< ^^"^^llj/bllAlU, Vj/GZ, VAgA, 

77 

Ca(A,A)<a||A||A||A||A, VA,AgA 


and of the linear form and/ 2 ,a: IKi,a||z' < {l-a)\\y,y^obs\\L^{rT) ||/ 2 ,a||A' < C(\\y^,obs\\L^(rT)- 

Moreover, with a G (0,1), we also obtain the coercivity of and of Cq,: precisely, we check 
that ar^a{y,y) > ^l|y||| for all y & Z, while Ca(X,X) > Q!||A||1 for all A G A. 

The result [H Prop 4.3.1] then implies the well-posedness of the mixed formulation 0 and 
the estimate p^ . □ 

The a-term in jCr,a is a stabilization term: it ensures a coercivity property of Er,a with respect 
to the variable A and automatically the well-posedness. In particular, there is no need to prove 
any inf-sup property for the application ba- 


Proposition 2.2 If the solution (y. A) G Z x L'^[Qt) of ^ verifies A G A, then the solutions of 
& and 0 coincide. 


Proof- The hypothesis of regularity and the relation (141 imply that the solution (y, A) G Z x 
LHQt) of § is also a solution of ( [l7| ). The result then follows from the uniqueness of the two 
formulations. □ 


Remark 2 Remark that the following unique continuation type property for 0 

^Ly = 0 m Qt) d^y = 0 on y = 0 in Qt, 

which is weaker than suffices to prove that the bilinear form (•,-)z in is a scalar product 
and therefore the well-posedness of We shall use specifically the observability inequality in 
the numerical Section Q to get an estimate o/ jjy - y/i||A 2 (Q^) from estimate of \\L{y - yh)\\L'^i^Q.^) 
and \\d,,{y - yh)\\L2{ST)- 
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2.2 Dual formulation of the extremal problem ([^ and remarks 

As discussed at length in [T3] for which we refer to detail, we may also associate to the extremal 
problem an equivalent problem involving only the variable A. Again, this is particularly 
interesting at the numerical level. This requires a strictly positive augmentation parameter r. 

For any r > 0, let us define the linear operator Vr from L'^{Qt) into LF‘{Qt) by 

VrX := Ay, VA e L^{Qt) 

where y ^ Z is the unique solution to 


ar{y,y) = h[y,X), yy&Z. 


(25) 


The assumption r > 0 is necessary here in order to guarantee the well-posedness of (25). 
for any r > 0, the form defines a norm equivalent to the norm on Z. 

We then have the following two results, proved in m Section 2.2] in the very close 
the exact controllability for 0 (see also Remark]^ below). 


Precisely, 
context of 


Lemma 2.2 For any r > 0, the operator Vr is a strongly elliptic, symmetric isomorphism from 
LP'{Qt) into LP‘{Qt)- 


Proposition 2.3 For any r > 0, let yo G Z he the unique solution of 

ar{yo,y) = l{y), yy & Z 

and let J** : L'^{Qt) —>■ K 6e the functional defined by 


Jr*W = -^ VrXXdxdt-b{yo,X). 
^ JJQt 


The following equality holds : 


sup inf £^( 2 /, A) = - inf J**(A) +£r(yo,0). 


This proposition reduces the search of y, solution of problem 0 , to the minimization of J**. The 


well-posedness is a consequence of the ellipticity of the operator Vr stated in Lemma 2.2 


Remark 3 We emphasize that the mixed formulation 0 has a structure very close to the one we 
get when we address - using the same approach - the null controllability of 0-' more precisely, the 
control of minimal L^(TT)-norm which drives to rest the initial data {yo,yi) € L‘^{Xl) x is 

given by v = lr,j, where (ip, A) G $ x LP'{Qt) solves the mixed formulation 


j a{(p,(p) + b{(p,X) = i{ip), 

1 %,A) = 0, yX€L^iQT), 

where ^ := Z, a and b are given by 0 and 0 respectively, while 1 is given by 

f : $ M, l{ip) = - J^pt{x, 0 )yo{x)dx + {g>{-,0),yi} 

and depends on the initial data (yo^yi) to be controlled. We refer to for the one dimensional 


case. 
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Remark 4 Reversing the order of ’’priority” between the constraint d^y — ,obs = 0 m Lf‘{TT) 
and Ly — f — 0 in L^{Qt), a possibility could be to minimize the functional y i—>■ \\Ly — /|li 2 (Q^) 
over y G Z subject to the constraint d^y — yu,obs = 0 in via the introduction of a Lagrange 

multiplier in However, the fact that the following inf-sup property: there exists i5 > 0 such 

that 

Ufj. d(^dt 


inf 


associated to the corresponding mixed formulation holds true is an open issue. On the other hand, 


> 6 


if a £-term is added as in (QRe ), this property is satisfied (we refer again to the book mu)- 


3 Recovering the source and the solution from a partial ob¬ 
servation: a mixed re-formulation of the problem 

Given a partial observation yv^obs of the solution on the subset Tt, we now consider the recon¬ 
struction of the full solution as well as the source term /. The situation is different with respect 
to the previous section, since without additional assumptions on /, the couple (y, /) is not unique. 
Consider the case of a source / supported in a set which is close to x {T}: from the finite 
propagation of the solution, the source / will not affect the solution y on Ft- 

Consequently, we assume that f{x,t) = a{t)fi{x) with a G C'^([0,r]) such that cr(0) ^ 0 and 
p G H~^{Ll) and then recall the following result of Yamamoto-Zhang (see [23 Theorem 2.1]): 

Theorem 3.1 Let us assume that the triplet {Vt-iT,Qt) satisfies the geometric optic condi¬ 
tion. Let y = y{p) G (^([O, Tj; Rq (fl)) n Tj; L^(r2)) be the weak solution of ^ with c := 1 

and (yo)2/i) = (OjO). Then, there exists a positive constant C such that 

C~^\\y\\H-^{n) < \\c{x) d,,y\\L2(rT) < C\\p\\h-^(Q), ^IJ^ e 

Therefore, assuming that the initial condition vanishes, this stability result implies the uniqueness 
of the source. 

We then consider the following extremal problem : 


[ inf J(y,/r) := -\\c{x){d,,y - y,,^obs)\\l 


1 




i'Py.u) 


subject to (y, p) GW 
where W is the space defined by 

W:= l^(y,p);yGC(lO,T];H((n))nC^(lO,T];L^(fl)), pGH-\Ll), 
Ly - ap = Oin Qt, y(-, 0) = yt(-, 0) = ol. 


(26) 


Note that, from Theorem 3.1 9^y S L‘^(Tt) for any {y,p) G W. We easily check that W is 
a Hilbert space endowed with the norm ||(y,/r)||w := I|c(x)i9,yy||i2(ri,). The extremal problem 
{'Py,u) posed : the functional J is continuous over W, is strictly convex and is such that 


J(y,p) -G - 1-00 as ||(y,/j,)||w —t oo. Moreover, in view of Theorem 3.1 the solution p is uniformly 
bounded in iJ“^-norm. 


As in the previous section, in order to avoid the resolution of the extremal problem (Vy^^i 


by an iterative process, we introduce a mixed formulation taking the equation Ly — ap = 0 as a 
constraint equality in L^{Qj’). In this respect, we define the following non empty vector space: 
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:= y g c([o,r];i7i(ii)) nci([o,r];L2(j7)), ^ g H-\n), 

Ly - an & L‘^{Qt), y(-, 0) = yt(-, 0) = o|. 

Then, as in Section we introduce the following hypothesis: 

Hypothesis 2 There exists a positive constant Cobs = C{Tt,T, ||c||^i(q)! II'^IIl“(0)) such that the 
following estimate holds: 

(^ 2 ) 


ff-i(n) < Cobs \\c{x)d^y\\L^^^^ + \\Ly- crnh^^Q^) ), V( 2 /, n) G Y. 


Again, in the case where the velocity c is constant, this inequality is a consequence of Theorem 
es to decompose any 

Lyi = ay in Qt, yi = 0 on Sy, 


3.1 it suffices to decompose any y with (y, y) gY as y = yi + 1/2 where yi and y 2 solve 

( 2 /i(-, 0 ),yi,t(-, 0 )) = (0,0), 


and 


Ly 2 = Ly-an in Qt, 2/2 = 0 on St, (y 2 (-, 0), y 2 ,t(-, 0)) = (0,0), 


respectively and then, applying Theorem 3.1 for yi and m Theorem 4.1] for y 2 , to write 

<l|c(x)(9,,yi||T2(r.r) < l|c(x)(9^y||T2(r,r) + ||c(x)a,,y2||T2(r,j,) 

< l|c(x)5i.y||T2(r.j,) + C{0,T)\\Ly - an\\mo,T-,L 2 {n)) 

for some (7(0, T) > 0. 

Then, within this hypothesis, for any y > 0 as in Section]^ we define on Y the bilinear form 

(( 2 /, m ), ( 2 /, m))v := // c^{x)di,ydoydadt + y {Ly - ay) {Ly - ajl) dxdt, ( 27 ) 
J j Y''j' »/«/ Q'j' 

for every {y,fJ.)i (y,M) G Y. We note the corresponding semi-norm ||(y,/i)||v := -\/((y,/.(), {y,n))Y- 
Again, in view of Q, d^y G L'^{Tt) for any (y,/i) G Y. We have the following result: 

Lemma 3.1 Under the hypotheses {TL 2 ), the space Y is a Hilbert space with the scalar product 
(•, ■)y defined by {21). 

Proof- From the inequality ( 7 ^ 2 ), the semi-norm || • ||y is indeed a norm. Let us check that 
Y is closed with respect to this norm. Let us consider a sequence {{yk, nk)}k>i G Y such that 
{yk,nk) —t ( 2 /, y) for the norm Y. Then, there exists / G Lf{Ll) such that Lyk — apik —t / in L'^{Qt) 
and g G L'^(Tt) such that d^yk g G ifiVT). Then, (772) implies that pk ^ y in 77“^(0). 
Consequently, Lyk = {Lyk — apk) + o’pk converges to f + ap in L'^{0,T; 77“^(0)). Consequently, 
yk can be considered as a sequence of solutions of the hyperbolic equation with zero initial data 
and second hand term converging in L^(0, T, i7“^(r2)). Therefore, by the continuous dependence 
of the solution, yk ^ y in C{[0,T]\L'^{O)) C (7^([0, Tj; i7“^(n)) where y is the solution of Ly = 
f + ap G L‘^{Qt) with (y(-, 0), yt(-, 0)) = (0,0) G V. But, again, from Theorem |3.l[ this implies 
that d^y G L'^{Vt) (and therefore g = d^y on Ft) and y enjoys the regularity (7([0, Tj; 77g (fl)) C 
C^{[0,T]-, L^{0.)) as the sum of two solutions (as for yi and y 2 above) in this space. Therefore, 

{y,p)GY. _ □ 

we introduce a Lagrange multiplier A G L^{Qt) and the following 


Proceeding as in Section 


2.1 


mixed formulation: find {{y,p),X) gY x L'^{Qt) solution of 

[ a{{y,p),{y,p)) Yb{{y,p),\) = l(y,p), 

\ &((2/,2i),A) = 0, 


V(y,/i) G Y 
WXGL^iQT), 


(28) 
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where 


a((?/,/i), (?/,/i)) := // C{x)d^yd^y dadt, 

JJtt 

b : Y X L‘^(^Qt) —^ I 

Z : yK, l{y,y,):= // c^{x) d^yy^^obadadt. 

JJtt 


:= // X{Ly - ay)dxdt, 
Qt 


Theorem 3.2 Under the hypothesis {H 2 ), the following holds : 


(i) The mixed formulation (23^ is well-posed. 

(a) The unique solution {{y,p),X) & Y x L‘^{Qt) is the saddle-point of the Lagrangian L : 
Y X L^{Qt) —>■ defined by 

^{{y,T),X) ■■= ^a{{y,fi),{y,p.))-\-b{{y,fi),X)-l{y,y). 


Moreover, the pair {y,p) solves the extremal problem {Vy^^}. 
(Hi) The following estimates hold : 


\\{y,H)\\Y = ||c(a:)(9^y||i2(r^) < \\c{x) y^^obsh^iTj 


and 


(29) 

(30) 


||A||l2(Qt,) < ^^/Cq^t + vWcix) y„ .o6s||L2(ry) 

for some constant Cq^t > 0. 

The proof is very close to the proof of Theorem |2.1| In particular, the inf-sup property can be 


obtained by taking, for any A £ L^(Qr), Mo = 0 and yo as in (12) leading to (j/OiMo) £ T so that 
the inf-sup constant 

b{{y,y),x) 


S := inf 


sup 


\eL^{QT) (y,y)(iY \\{y^ y)\\Y\\X\\L2{Qj.) 
is bounded by above by {Cq^t + r])~^^’^. 

Then, as a consequence of the estimate [Hi I, the source term y is uniformly bounded in the 
iJ~^-norm: 


Corollary 3.1 Assuming (Hi), the solution p of (28) satisfies\\yL\\fi-i^ri) Si Cobs\\c{x) y,^^obs\\L'^{rT 

Again, it is very convenient to ’’augment” the Lagrangian (see [T^) and consider instead the 
Lagrangian Cr defined for any r > 0 by 

Criiy, y),X) := ^ar((y, p), {y, p)) Y b{y, X) - l{y, p), 

aT{{y,p), (y,M)) := a{{y,p), {y,p)) + r\\Ly - apWl^t^Q.^). 

Since a{{y, p), {y, p)) = ar{{y, p), {y, p)) in W, the Lagrangian £ and £r share the same saddle- 
point. The positive number r is an augmentation parameter. Moreover, proceeding as in Section 
we may also associate to the saddle-point problem sup;^g^ 2 (Q,j,) inf(y^^)gy Cr{{y, p), A) a dual 


2.2 


problem, which again reduces the search of the couple {y,p), solution of problem (Vy^y), to the 
minimization of a elliptic functional with respect to A. 
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Proposition 3.1 For any r > 0, let (yo,Aio) &Y be the unique solution of 

ar((2/o,Ato),(|/,M)) = V(y,7i) e Y 

and let Vr be the operator from L‘^{Qt) into L‘^{Qt) defined by := Ly — cr/i where {y,p,) G Y 
is the unique solution to 


'i{y,n)GY. (31) 

The operator Vr is strongly elliptic and symmetric. Moreover, the following equality holds 

sup inf Cr{{y,p,),X) = - inf + Cr{{yo, y.o),0) 

\eL^{QT)0,u)&y 

where Jf* : L^{Qt) —>■ M is the functional defined by 


J:*{X) = 1 jj {VrX)Xdxdt-biiyo,fio),X)- 


Proof- From the definition of Or, we easily get that ||’Pr-A|ji 2 (Q^) < ||A|ji 2 (Q^) and the 

continuity of Vr- Next, consider any A' € L^{Qt) and denote by {y',yl) the corresponding unique 
solution of (31| so that VrX’ := Ly' — ap'. Relation (31) with {y,Jl) = {y', p') then implies that 

(32) 


[[ {VrX')Xdxdt = ar{{y,p),{y',p')) 
JJqt 


'Q 

and therefore the symmetry and positivity of Vr- The last relation with X' = X and the observability 
estimate ("^2 I imply that Vr is also positive definite. 

Finally, let us check the strong ellipticity of Vr, equivalently that the bilinear functional 

(A, A') —>■ // {VrX)}f dxdt is L^((5T)-elliptic. Thus we want to show that 
JJqt 


[[ {VrX)Xdxdt>C\\X\\l.^Q^y ^XgL^Qt) 
J J Qx 


(33) 


for some positive constant C. Suppose that (33) does not hold; there exists then a sequence 
{A„}„>o of L'^{Qt) such that 


ll-^nllL2(QT) = Vn>0, lim // [VrXrfjX^dx dt = Q. 


Let us denote by {y„, Pn) the solution of (311 corresponding to A„. From (32), we then obtain that 
lim (r||Lj/„ - cTPn\\l2(Q.r) + l|c(a:) duyn\\l 2 (rr.)) = 0- (34) 


From (31) with {y, p) = {yn, Pn) and A = A„, we have 


[r{Lyn - apn) - Xn]{Ly - ap)dxdt+ c{x)d,,ynduyda dt = Q, \/{y,p)GY. (35) 

Q'X' J J T''J' 

We define the sequence {y„}n>o as follows : 

Ly„ = r {Lyn - apn) - A„, in Qt, 

Vn = 0 ) on St, 

yn(')O) = y„_t(-,0) = 0, in O, 

so that, for all n, is the solution of a hyperbolic equation with zero initial data and source 
term r (L?/„ - apn) - A„ in L'^{Qt). Using again [3T1 Theorem 4.1], we get \\dyyr,\\L‘^{TT) < 
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\/C^\\r{Lyn - a^n) - 
we get 


so that (y„, 0) S Y. Then, using (351 with {y, y) = (y„, 0) S F 


\\r{Lyn - cr/x„) - KWh^iQ^.) < \/C'a,T||c{a;)5j,j/„||z,2(ri,). 


Then, from (34), we conclude that lim„_>.+oo ll^n||L 2 (Q,j,) = 0 leading to a contradiction and to the 
strong ellipticity of the operator Vr- The rest of the proof is standard. □ 

To end this section, it is worth to mention that a result similar to Theorem 3.1 is given in 
|261 Theorem 1] with a control of y in L^{Lt). However, this requires an observation of d^y in 
T; L^(r)) (equivalently of dtd^y in L^(rj’)) which may be too strong in real life applications. 


4 Numerical Analysis of the mixed formulations 

4.1 Numerical approximation of the mixed formulation (|^ 

We now proceed to the numerical analysis of the mixed formulation ([^, assuming r > 0. We follow 
[13], to which we refer for the details. 

Let Zh and A/j be two finite dimensional spaces parametrized by the variable h such that 
Zh C Z,Ah C L'^{Qt) for every h > 0. Then, we can introduce the following approximated 
problems: find {yh,^h) & Zh x Ah solution of 


ariyh,yh) + KyhAh) 

b{yh,\) 


Kyn)^ 

0 , 


'^yh G 
yXh G Ah- 


(36) 


The well-posedness of this mixed formulation is again a consequence of two properties. The first 
one is the coercivity of the bilinear form on the subset Afh{b) = {yh G Zh', b{yh, Xh) = 0 VA/i S 
A/i}. Actually, from the relation ar{y,y) > (r/? 7 )||y ||2 for y ^ Z, the form ar is coercive on the 
full space Z, and so a fortiori on Nh{b) C Zh <Z Z. 

The second property is a discrete inf-sup condition: precisely, for any h, 


Sh := inf sup 


b(yh,Xh) 


\\Xh\\L^(QT)\\yh\\z 


> 0 . 


(37) 


Let us assume that this condition holds, so that for any fixed h > 0, there exists a unique 
couple {yh,Xh) solution of (36). We then have the following estimate. 


Proposition 4.1 Let h > 0. Let {y, X) and {yh, Xh) be the solution of & and of (g6| ) respectively. 
Let Sh be the discrete inf-sup constant defined by (Sly. Then, 


\y-yh\\z < 2 i-H 




d{y,Zh) -L -—d{X,Ah), 

Vv 


IIA - A.|U2(q,) <(^2 + 1%, Zh) + A,) 

where d{X, Ah) := infA,.GAfe |1 A - XhWL'^^Q,^) 


(38) 

(39) 


d{y,Zh) 


:= inf \\y-yh\\z 

Vh&Zh 


= inf 



\ 1/2 

9t^yh\\'L2(rT) + v\\L{y — i//i)|Il2(q,j,) j 


( 40 ) 
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Proof- From the classical theory of approximation of saddle point problems (see [H Theorem 
5.2.2]) we have that 


\y-yh\\z < 


2\\a 


+ 


and 


r\\(ZxZ)' 
-T 

ao 

IHjZxL^iQT))' 

ao 


2||a, 


r\\(zxzy\o\\{zxL^QT)y 


anS, 


d{y,Zh) 


0 


d{\,Ah) 


l|A - < 


2||a. 


r\\{ZxZy 


■r||(ZxZ)'l|0|l(ZxL2(Qy))' 


an (5/1 


51 


d{y, Zh) 


3||ar|| 2 I1 ^II(ZxL2(Qt))' 

1 

anSh 


d(A,A„). 


Since, \\ar\\(zxzy = ao = 1; 


zxL2(Qt))' = the result follows. 


Remark 5 For r = 0, the discrete mixed formulation (36) is not well-posed over Zh x Ah because 


the form ar=o is not coercive over the discrete kernel ofb: the equality b{yh, Xh) — 0 for all Xh € Ah 
does not imply in general that Lyh vanishes. Therefore, the term r||L?//i||| 2 (Q,j,), which appears in 
the Lagrangian Cr, may be understood as a stabilization term: for any h > 0, it ensures the uniform 
coercivity of the form and vanishes at the limit in h. We also emphasize that this term is not 
a regularization term as it does not add any regularity on the solution yh- 

Let nh = dmi{Zh),mh = dim(A/i) the dimension of the space Zh and Ah respectively. Let the 
real matrices Ar^t € K"'*’”'*, Bh G Jh G and Lh G K"'* be defined by 


ar{yh,yh) = (Ar,/i{?//t}, {yft,})R’*fe,R"h 
b{yhiXh) = {A/i})R™h,R"*fe 

XhXh dx dt — 


Qt 

. l{yh) = {Lh,{yh})w 


'dyh,yh G Zh, 
'dyh G Zh, Xh G Ah, 

'dXh, Xh G Ah, 

'dyh G Zh, 


(41) 


where {yh} G K"** denotes the vector associated to yh and (•, •)R"fe,H"h the usual scalar product over 


‘. With these notations, the problem (36) reads as follows: find {yh} G M"'* and {A/j} G K'"'* 
such that 


A. 


r,h 


0 


(42) 


{yh} 

^ R"h+™fi>"h+"*h V {A/l} / V " / 

The matrix Ar^h as well as the mass matrix Jh are symmetric and positive definite for any h > 0 


Lh 

0 


and any r > 0. On the other hand, the full matrix of order mh -\- nh in (42) is symmetric but not 
necessarily positive definite. 


We recall (see [H Theorem 3.2.1]) that the inf-sup property (37) is equivalent to the injective 
character of the matrice Bj( of size nh x mh, that is Ker{Bj() =0. If a necessary condition is 
given by mh < Uh, this property strongly depends on the structure of the spaces Zh and Ah- We 
will discuss numerically this property in Remark for a specific choice of approximation. 

4.1.1 (7^-1111116 elements and order of convergence for A^ = 1 

The finite dimensional and conformal space Zh must be chosen such that Lyh belongs to L‘^{Qt) 
for any yh G Zh- This is guaranteed, for instance, as soon as yh possesses second-order derivatives 
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in L‘^{Qt)- As in [O] . we consider a conformal approximation based on functions continuously 
differentiable with respect to both variables x and t. 

We introduce a triangulation Th such that Qx = UxeTh^ assume that {Th}h>o is a 

regular family. We note h := max{diam(iir), if G Th} where diam(it') denotes the diameter of K. 
Then, we introduce the space Zh as follows : 

Zh:={zh€C^(^):zhlK€P(K) VK € Th, Zh = 0 on Ex}, (43) 

where P(i^) denotes an appropriate space of functions in x and t. In this work, we consider two 

choices, in the one-dimensional setting (for which C M and Qx C K^): 

(i) The Bogner-Fox-Schmit (BFS for short) C^-element defined for rectangles. It involves 16 

degrees of freedom, namely the values of yh,yh,x,yh,t,yh,xt on the four vertices of each rect¬ 
angle K. Therefore ¥{K) = P 3 ® P 3 ^ where Pr,^ is by definition the space of polynomial 

functions of order r in the variable We refer to [TUI ch. II, sec. 9, p. 94]. 

(ii) The reduced Hsieh-Clough-Tocher (HCT for short) C^-element defined for triangles. This is 
a so-called composite finite element and involves 9 degrees of freedom, namely, the values of 
yh,yh,x,yh,t on the three vertices of each triangle K. We refer to [ini ch. VII, sec. 46, p. 
285] and to jU] [53] where the implementation is discussed. 

We also define the finite dimensional space 

Ah := {A,, e c°(^), \h\K e Q{K) VK G %} 

where Q{K) denotes the space of affine functions both in x and t on the element K. 

For any ft, > 0, we have Zh C Z and Ah C L^(Qt)- 
We then have the following result: 

Proposition 4.2 (BFS element for iV = 1 - Rate of convergence for the norm Z) Let 
h > 0, let k < 2 be a nonnegative integer. Let (y. A) and {yh,^h) be the solution of & and 
respectively. If the solution (y. A) belongs to H^^^{Qx) x H^{Qx), then there exist two positives 
constants 

Ki = Arj(||?/||^fc+2(Q,j,), ||c||c.i(^), ||d||i=o(Q,j,)), i G {1, 2}, 
independent of ft, such that 

\\y-yk\U<Ki(l + ^^ + Ly\ (44) 

IIA - A4|U.,a,, < A'2 ((i + I + (45) 

Proof - From [TU] ch. Ill, sec. 17], for any A € H^{Qx), k <2, there exists Ci = Ci(||A||^fc(Q,j,)) 
such that 

||A-nA,.r.(A)||L=(Q.) Vft>0 (46) 

where IIa^^^t-,^ designates the interpolant operator from L‘^{Qx) to Ah associated to the regular 
mesh Th- Similarly, for any y G H^~^^{Qx), there exists C 2 = C2{\\y\\H>‘+^{QT)) every 

ft > 0 we have 


\\y-^Zh,rhiy)\\L^iQT) ^ C2h'"~^'^, \\y-iizh,rh(.y)\\H^{QT) ( 47 ) 

Then, using that the linear operator P : H^T(^Qrp') —j. LflEx) defined by Py := d,,y is continuous, 
there exists a positive constant Cp such that 


\\c{x){d„y - d,,(nz,^,Thiy)))\\L^{rT) < \\c\\L°°{rT)Clp\\y - nz,^,Th{y)\\m/^{QT)- 


(48) 
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We then observe that 

\\Ly- L{Uz^,rh{y))\\L^QT) < ^(llc|lc’i(^)’ IMIlL“(QT))lly “ ^Zh,n{y)\\H^QT)^ 
for some positive constant K; ( |40| then leads to 

d{y,Zh) < {\\d^y - d^{Uz,,,TH{y))\\LHrT) + v\\Ly - L{nzh,n{y))\\L^QT)] 

and then from Proposition |4.1[ we get that 


(49) 


(50) 


\y-yh\\z <21 




Similarly, 


IIA - < (2 + ^)yil4U,^,C}+vK-C,k'- + 

From the last two estimates, we obtain the conclusion of the proposition. □ 

It remains now to deduce the convergence of the approximated solution for the L‘^[Qt) 
norm. This is done using the observability estimate Precisely, we write that {y — yh) solves 
the hyperbolic equation 

L{y - yh) = -Lyh in Qt 

{iy-yh),iy-yh)t)i-, o) e y 

y-yh = 0 on St. 

The continuous dependence of the solution with respect to the right hand side and the initial data 
leads 

\\y-yhWhiQ^) - C!n,T{\\{{y - yh)i-,0), (y - yh)t{-,0))\\v + WLyhWhiQ^))- 
Combined with @ applied for y — yn, this leads to 

\\y ~ yh\\'L'i(Q^) < Cn,T{7 + Cobs){\\y — 2/?i|li2(r,j,) + \\^yh\\'L'i(QT)) 

from which we deduce, in view of the definition of the norm Y, that 


\\y - yhWh^QT) < C'n.T(l + C'of,s)max(l,—)||y-?//,||z. 


(51) 


Eventually, by coupling (51) and Proposition 4.2 we obtain the following result 


Theorem 4.1 (BFS element for N = l - Rate of convergence for the norm L'^{Qt)) 


Assume that the hypothesis 


holds. Let h > 0, let k < 2 be a positive integer. Let [y, A) and 


(.yh,^h) be the solution of and (36) respectively. If the solution {y,X) belongs to H^'^^{Qt) x 
H^{Qt), then there exists two positives constant K = R(||j/|lfffc+ 2 (Qj,), ||c||(;;,i(^), ||d|!L~(Qi,), Cq.t, C'obs), 
independent of h, such that 


lls- ^)(i + + <)a‘. 


(52) 


Remark 6 (About the choice of the augmentation parameter r) Estimate (52) is not fully 
satisfactory as it depends on the constant Sh ■ In view of the complexity of both the constraint Ly = 0 
and of the structure of the space Zh, the theoretical estimation of the constant 6h with respect to 
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h is a dijficult issue. However, as discussed at length in m Section 2.1], Sh can be evaluated nu¬ 
merically for any h, as the solution of the following generalized eigenvalue problem (taking ij = r, 
so that ar{y,y) is exactly ||j/||U.' 


Sh = inf<j : BnAflBliXh} = 6Jh{Xh}, 


V{A4gM'"'*\{0} 


(53) 


where the matrices Ar^h> and Jh are defined in m- 

The computation of Sh given by the problem (53) with respect to h and to r = rj has been 
performed in Section ].2] for the following data : T = 2, Qt = (0,1) x (0,T), T = 2, 
= {1} X (0,T). There, it is observed, that for both the BFS and the HCT finite element on 
regular meshes (non necessarily uniform), the constant Sh behaves like 

Sh « Cr—}= as /i —>■ O’*' (54) 

\ r 


with Cr > 0, a uniformly bounded constant w.r.t. h. For the BFS finite element, Table^ reports 
some numerical values of Sh while depicts the evolution of ^/rSh w.r.t h for r = 1,10“^, h 

and r = hf. 


h 

1.41 X 10-1 

7.01 X 10-2 

3.53 X 10-2 

1.76 X 10-2 

8.83 X 10-3 

r = 10^ 

4.71 X 10-2 

2.53 X 10-2 

1.23 X 10-2 

6.06 X 10-3 

3.01 X 10-3 

r = 1 

4.17 X 10-1 

2.53 X 10-1 

1.23 X 10-1 

6.05 X 10-3 

3.01 X 10-2 

r = 10“^ 

1.474 

1.427 

1.178 

5.99 X 10-1 

3.01 X 10-1 

r = h 

1.301 

1.086 

7.65 X 10-1 

5.37 X 10-1 

3.79 X 10-1 

II 

to 

1.665 

1.483 

1.485 

1.489 

1.497 


Table 1: Sh w.r.t. r and h - T = 2. for the BFS element. 



Figure 1: BFS finite element - Evolution of y/rSh^r with respect to h (see Table for r = 1 (□), 
r = 10“^ (o), r = h (k) and r = h^ (<). 


Reporting this behavior in the estimate (52), the error \\y — yh\\L^{QT) behaves like (taking rj = r) 
\\y - VhhRQT) < K max(l, + \ + ■ 
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The right hand side is minimal for r of the order one leading to \\y — yh\\L^{Qj,) ^ Kh^~^. This 
estimate is very likely not optimal but shows the strong convergence of the approximation for regular 
enough solution. In particular, the estimate (4-8) of the boundary term by the H^{Qx)-norm may 
result in a loss of precision. 

On the other hand, the same argument for the variable Xh indicates, using (45) that the estimate 
of ||A - is 


l|A - A/j||i2(Q^) < K2'^(^ 




— \h‘" 


The optimal value for the augmentation parameter is now r = hf leading to ||A — A;i||j;, 2 (Q^) < 
K 2 h^~^. Remark that for r = , the discrete-inf sup constant 5h remains uniformly bounded by 

below w.r.t. h. As observed in m Section 4-] (see also the numerical experiments of the present 
paper), if the influence of the parameter r > Q is not really important (for h small enough), the 
choice r = h^ offers the best rate of convergence. 

As discussed and used in Section 4-3], when y^ is approximated with the BFS and HCT finite 
element, the quantity h\\Lyh\\L^(^Q.j.) is asymptotically equivalent w.r.t. h to \\Lyh\\L^{p^T-,H-'^{n))- 
Therefore, taking r = h^ in the augmented mixed formulation amounts to relax the constraint 
Ly = 0 in L^{Qt) by the weaker one Ly = 0 m T, H~^{Ll)) enough in practice to approximate 
weak solutions of 0 - 


4.2 Numerical approximation of the stabilized mixed formulation ( |17[ ) 

We address the numerical approximation of the stabilized mixed formulation ( |17[ ) with a € (0,1), 
assuming again that r > 0. Let h he a real parameter. Let Z/i and be two closed finite 
dimensional spaces such that 


Z/, c Z, AhC A, V/i > 0. 


The problem (17) becomes : find {yh,Xh) G Zh x Ah solution of 


ar,a{yh,yh)+ bai>^h,yh) = hAVh)^ Ah^Zh 

bct{Xh, yh) r-aifXh, Xh) — ^2,a(A/j), VA/j G Ah- 


(55) 


In view of the properties of the forms Ca, and ^ 2 ,q) this formulation is well-posed. 

Proceeding as in the proof of [H Theorem 5.5.2], we first check that the following estimate 
holds. 


Lemma 4.1 Let {y,X) GY x A be the solution of (11) and (yh,Xh) G Zh x Ah be the solution of 
(55). Then we have. 


— UhW'z + ^oi\\X — Xh\\\ < 


'zxzy ll^“ll(zxA)' 0 


9 


Y 


+ A )_inf hh-yWz 
2 J 


ll^“ll(ZxA)' , , ■ r iiT mi2 

- + a + - _mL IIA,, - Alii 

^ ^/AfcGAfe 


(56) 


with ||ar,a||(zxz)' < max(l — a,rj ^r), || 6 a||( 2 xA)' < {Co,,t + otfr] Parameter 0 is defined in 
Proposition \2.1\ 

Concerning the space Ah, since LXh should belong to L'^{0,T,H~^{D)), a natural choice is 

A;, = {AeZ^;A(-,0) = Ad-,0) = 0}. (57) 


where Zh C Z is defined by (43). Then, using Lemma [4.1| and the estimate (50), we obtain the 
following result. 
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Proposition 4.3 (BPS element for N = 1 - Rate of convergence for the norm Z x A) 
Let h > 0, let k <2 be a positive integer and a € (0,1). Let (y, A) and {yh, Xh) be the solution of 
(l ?{ ) and (55) respeetively. If (y, A) belongs to x then there exists a positive 

constant K = i4r(||y||^fe+2(Q^), a, r, y) independent of h, such that 

\\y-yh\\z + \\X-Xh\\A<Kh^ (58) 

Sketch of the proof- The estimate of ||y — Hzh.,Thiy)\\z for any y G Z in term of 0{h^)^ is 
detailed in the proof of Proposition 4.2 precisely, we refer to (50). Similarly, we write that, for 
any A € A, 


inf ||A-A^||a< ||A-n^^_^JA)|lA 




< ( ||L(A - n^,^_ 7 -,^(A))||| 2 (^-i(Q)) -I- ||c(j;)(A - n^,^_ 7 -,^(A))||| 2 (r.j 


1/2 


< (^||L(A)-L(%^_^^(A))||i2(Q,) + ||c||i,.(r.)l|A-%^_^^ 

< Ki\\c\\c^^)’ II^^IU“(Qt))I|A - nA;.,rfc(A)llH2(QT) 


(A)||^ 1 / 2 (q^- 


1/2 


(59) 


for which the 0{h^) estimate follows, in view of the definition of A/,. 

In particular, arguing as in the previous section, using Proposition |4.3[ the observability esti¬ 
mate for the variable y and the estimate ( [T^ for the variable A, we get the following desired 
global estimate: 

Theorem 4.2 (BPS element for A^ = 1 - Rate of convergence for the L'^{Qt) norm) 
Assume that the hypothesis & holds. Let h > 0, let an integer k <2. Let (y. A) and {yh, \h) be the 
solution of ( f/7| ) and (55) respectively. If the solution (y. A) belongs to H^^^{Qt) x H^^‘^{Qt) , then 
there exist two positive constants K,, = Ki{\\y\\Hi=+ 2 (^Q.^), ||A||/i-fe+ 2 (Q,j,), \\c\\ci(^), ||d|| a=o(q^), a, r, y), 
/ = 1,2 independent of h such that 

h^ 


\\y - yhWh'^iQT) < ~ Aa||l2(Qt) < K2h'". 


(60) 


We emphasize that, within the stabilized formulation, these estimates do not depend of any discrete 
inf-sup constant. In particular, the positive augmentation parameter r can be chosen arbitrarily. 


4.3 Numerical approximation of the mixed formulation (28) 


We now consider the numerical analysis of the mixed formulation ( |28[ ) where both the solution 
y and the spatial part y of the source are unknown. We take a strictly positive augmentation 
parameter r. 

Let Yh and Ah be two finite dimensional spaces parametrized by the variable h such that 
Yh C Y,Ah C L'^{Qt), and dim(yA) > dim(A^) for every h > 0. Then, we can introduce the 
following approximated problems: find {{yh, y-h), Xh) GYhX Ah solution of 


I ar{{yh,yh),{yh,yh))+ b{{yh,yh),Xh) = i{yh), '^{yh,Th)^Yh 

I b{{yh, yh),Xh) = 0 , VAa G A;,. 

Again, for any r > 0, the coercivity of the form Or holds true on the whole space Y: 

ar{{y,y),{y,y)) > ^{y,y)\\Y, V(y,/i) g Y 

V 


( 61 ) 
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so a fortiori, uniformly on the subspace Yh- Therefore, assuming that the discrete inf-sup constant 


defined by 


Sh ■■= inf 


sup 


K[yh,t^h)Ah) 


\\^h\\L‘^(QT}\\{yh, yh)\\ 


Y 


(62) 


is, for any h 
holds true: 


> 0, strictly positive, there exists a unique solution {yh,^ih) of (611 and the following 


Proposition 4.4 Let h> 0. Let and {{yh, fJ-h), Xh) be the solution of \28) and of {Sf 

respectively. Let Sh be the discrete inf-sup constant defined by Then, 


\\{y,T) - {yh,yh)\\Y < 2 i -f 


l|A - Xh\\L^{Q.r) ^ 


1 




1 


d{{y,p),Yh) + ■^d{X,Ah), 


Yn^h J Sh 

where d{X,Ah) is as in Proposition^^ and 


—d{y,Zh) 


Vvdi 


-d{X,Ah) 


(63) 

(64) 


d{{y,y),Yh) ■■= inf \\{y,y)- {yh,yh)\\Y 

{yh,^h)^yh 


i - yh)\\ l^tt) + y\\d-(y - yh) - - yh)\\L^Qj.) 


1/2 


(65) 


Proof- The proof is similar to the proof of Proposition 4.1 using again that ||ar||(yxT)' ^ 1 ami 
li&li(YxA)' < 

The finite dimensional problem (61) reads as follows: find {yh,yh} € R”'* and {A/j} G K'"'* 
such that 


A. 


r.h 


Bh 


K 

0 




{yh,yh} 

{A4 




Lh 

0 


( 66 ) 


with matrices similarly defined as in the previous section. 

Since the variable /r is a function of x S SI only, we first introduce a triangulation Ta^, such 
that SI = and we assume that {7Aa;}Ax>o is a regular family. We denote Ax := 

max{diam(Ar), K G Tax}- We then introduce the subspace Max of H~^{Ll) defined by 


Max — {t^Ax G1 C^{fl) : y,Ax\K G Qx{K), Vif G Tax}, 


(67) 


where Qx{K) denotes the space of affine functions in the variable x in the element K. 

Similarly, we define a subdivision Taj of [0,T] such that [0,r] = UggTAtQ and denote At := 
max{diam((3), (5 € Ta*}- Then, we consider the triangulation defined by {Th} ■= {7Ax}Aa:>o ® 
{7At}At>o such that 

Qt = [J {K X Q}. 

Again, we denote by h := max{diam(Ar x Q),K G Tax,Q G Tai}- The triangulation {Th}h is a 
regular family for Qt as soon as {7Aa;}Ax>o is a regular family for SI. 

As in the previous section, we now go on in the one dimensional case in space. 


4.3.1 C^-finite elements and order of convergence for N = 1 

We only discuss the BPS finite element and introduce the space Zh as follows : 

Zh = {zh C C^{Qt) : Zh\M G P(M), VM GTh, Zh = 0 on Ft, Zh = Zhu = 0 on SI x {0}}, (68) 
with P(M) := P 3 , 3 :(Q) x for all {Q,K) G (7Ax,^Ai)- Finally, we define the space Yh by 

Fk •— {y/i ~ f^Ax) • € Zh^l^Ax € -^Aa;} 
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so that, for each value of t G [0,T], the variables Zh{-,t) and fi^x share the same triangulation with 
respect to the variable x. We check that Yh C Y for each h > 0. In the sequel, for simplicity, we 
use the notation for /tax- 

We then have the following result: 


Proposition 4.5 (BPS element for N = 1 - Rate of convergence for the norm Y) Let 
h > 0, let k,q<2 be two nonnegative integers. Let and ({yhi k-h)iXh) he the solution of 

(28) and (61) respectively. If the solution ((y,/r),A) belongs to x x 

then there exist two positive constants 


Kt — ||/r||/i-<!(a), ||o'||l=x=([o,t])), * G {1,2}, 

independent of h, such that 


\iy,h') - iyh,yh)\\Y < -ffi (i + 

l|A - Xh\\L^(Q.^) < K 2 


1 


1 


h'^+KAl + 




*) 

1 

5h ^/rjdh 




(Aa 


1 


h'^+K^i 1 + 




ML 


(69) 

(70) 


Proof - We proceed as in the proof of Proposition |4.2| We write that 


\\L{y - IlYh,Thiy)) - - ^ZH,rh{T))\\L^(QT) 

< \\L{y - HzH.Th{y))\\L-^{QT) + lklU”([o,r])||M - ^M^:,,T^AT)\\L'^(n) 

< K\\y - Hzh,Th{y)\\H‘^{QT) + lkllL“([o.T])l!M - IImax.Ta, (li)IU2(n) 


where K is from (49). 

We use that, if G H'^{Ll), q < 2, there exists a positive constant C 3 = C'3(||/i||ij9(n)) such 
that ll/i — A C' 3 (Ax)‘^ for any Aa; > 0. Consequently, from (47), we have the 

estimate 


WHy - ny,.r.(y)) - ^(/a - nz^nipMLHQ^) < KC^h^ + C^iAx^ < {KC^ + C^Kh^^ + (Aa;)«). 
This leads to the following estimates : 

1/2 

) 

< V2^(||c||i^(r,)C2 + V2^\\ah^([o,T])C3{Axy 


d{{y,y),Yi,)< + 


Using (63), we then get 

\\iy,T)-{yh,p.h)h 


<2V2(l + y{\\c\\l^^^^^Cl + qK^)C2h>^ + yty\\ah^ao,T]))C3{Axr 


V Vy^h 
-^c,h\ 


Similarly, using (ro4|, we get 


l|A - Xh\\L^(QT.) 

<V2(^2 + 

3 


y/nShJ Sh 


+ yK^)C2h'^ + ^\\ah^ao,T])C3iAxy 
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The proposition then follows from the last two estimates. □ 

It remains now to deduce the convergence of the approximation yh for a global norm, typically 
We write that y — yh solves the problem 

L{y - Vh) = - yh) + L{y - Vh) - cr{y - yh) in Qt 

iiy-yh)Ay-yh)t)i-,o) = {o,o) infi 

y-yh = 0 on St 

leading to 

\\y - yhWL^iQr) < C'n.rf lkllL“(o.T)||(Ai - yh)\\H-^(n) + \\L{y - yh) - cr{y - M/i)IIl2(Qt) )• 


On the other hand, proceeding as before, assuming (7^2), we have 


\\y yh\\H-^{Q) — ^obs\^\c{x){d,^y 9L,yh)\\h^(j'j,) + \\L{y yh) cr^y M/t)llL 2 (Qr) 

leading to 

\\y - VhWh'^iQT) < C'a.T^lk||Loo([o_T])\/^^l|c(a:)(9^(?/ - yh)\\L^(VT) 

+ (1 + ||cr||L=°([ 0 .T])\/Cobs)|k( 2 / - Vh) - cr{y - yh)\\L^Q._ 

< Cq,t{I + lk||Loo([0,T])\/Cobs) (^\\cix)d^{y - yh)\\L^{rT) 

+ \\L{y - yh) - cr{y - yh)\\mQT)^ 

< k3C7n,T(l + lk||Too([o,T]) VC'o&s)niax(l, -^)\\iy -yh,y- yh)\\Y- 


Therefore, in view of Proposition 4.5 we get the following a priori estimate for the L^-norm : 

Theorem 4.3 (BPS element for = 1 - Rate of convergence for the L^((5t)-norm) 
Let h > 0, let k,q <2 be two nonnegative integers. Let ((j/,/r), A) and {{yh, yh), ^h) be the solution 
of (28) and {61} respectively. If {{y,y),X) belongs to {H^'^'^{Qt) x H'^{LL)) x H^{Qt), then there 
exists a positive constant 

K = Ar(||?/||^fc+2(Q,j,), |kllH''(a)i Ikllci(^); lkllL“(Qr)); 
independent of h, such that 

1 


\y - yhWL^iQr) < + lkllL“([0.T])\/Co&s) max(l, —) 

y/V 


1 + 


VySh ydFj 


^ 1 + 


Vv^h 


(Aa 


Remark 7 ITe have presented along Section 2] error estimates in the one-dimensional case in 
the case where the space Zh (and Ah for the present section) is based on the BFS finite element. 
Similar results may be obtained within the HCT finite element. Precisely, we refer to HE ch. VII, 
sec. 48, p. 295] where the following interpolation estimates for the HCT element are developed: 

\\y-Hz„,TM\\L-iQP < Ch^^\ \\y -^z„,rAy)\\HVQr) < Ch\ Vy^H^Qr) 

for k G {0,1} and some constant C > 0. 
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5 Numerical experiments 


We now report and discuss some numerical experiments corresponding to mixed formulations (361, 
(55) and (61 1 for TV = 1 and N = 2 respectively. 


5.1 Reconstruction of the solution - One dimensional case (N = 1) 

We take Lt = (0,1) and T = {1} C dU, T = 2, c = 1 and d = 0. We check that for these 
data, the inequality holds true. Moreover, in order to check the convergence of the numerical 
approximations, we consider explicit solutions of ([^. We first define the following initial condition 
in Hq{LI) X (see [9]): 


(EXl) j/o(a;) = 1 - |2x - 1|, yi{x) = l(i/ 3 , 2 / 3 )(a;), 
and / = 0. The corresponding solution of ([^ is given by 


e(o,i) 


y{x,t) = Uk cos{kTrt) + ^ sm{kTrt) )-^2 sin(/c7rx) 
k >0 ^ 


(71) 


with 


4^2 

Ofc = sin(7rfc/2). 


In particular WvWl^Qt) = 


bk = —(cos(7rfc/3) — cos(27rfc/3)), 
irk 


k>0. 


k >0 


+ Ji. 

k ' (kn)- 


) « 0.58663. The corresponding normal derivative 


dt^ylrT = yxi^U) is depicted in Figure]^ This example is rather stiff: the normal derivative is 
in L^(0,T) but is discontinuous (in view of the regularity of the initial condition). We compute 


Wdi^yh^rr) = \J‘^Y.k>oii^'^y^k + K) ~ 8.33298. 


We recall that the direct method amounts to solve, for any h, the linear system (42). We use 


exact integration methods developed in m for the evaluation of the coefficients of the matrices. 


Moreover, the linear system (42) is solved using the LU decomposition method. 


We first consider the BFS finite element with uniform triangulation (each element K of the 
triangulation 7)i is a rectangle of lengths Aa; and At so that h = (Ax)^ + (At)^). Table [^collects 

some numerical values with respect to h for r = h'^ and for Ax = At. We observe the following 
behavior with respect to h: 




Wy-yhWh^Q 

r = h'^ : WvWl^Qt) 

\\^h\\L2{QT) = 0 (^^'“), 


\\d,^iy - yh)\\LHrT) 


\\duy\\L^{rT) 




The evolution of the norm \\y — yh\\L2{QT) suggests that the unique solution y is correctly recon¬ 
structed from the observation. Moreover, in agreement with Remark since y,^^obs is by construc¬ 
tion the restriction to {1} x (0, T) of the normal derivative of a solution of 0 , we check, that the 
sequence A/i, approximation of A, vanishes as —>■ 0. Eventually, in view of the Remark]^ we have 
A\Lyh\\L'^{QT) = h\\Lyh\\L'^(Q.^) « ||iy?i||L2(0,T.H-i(0,i)) ~ so that yt approximates 

correctly the unique corresponding weak solution of ([^. 


We also check that the minimization of the functional J** introduced in Proposition 2.3 leads 


exactly to the same result approximation: we recall that the minimization of the functional J** 
corresponds to the resolution of the associated mixed formulation by an iterative Uzawa type 
method. The minimization is done using a conjugate gradient algorithm (we refer to |131 Section 
2.2] for the algorithm). Each iteration amounts to solve a linear system involving the matrix 
which is sparse, symmetric and positive definite. The Cholesky method is used. The performance 
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of the algorithm depends on the condition number of the operator Vr- precisely, it is known that 
(see for instance m), 


II A” 


A||l2(Q,j,) < 2s/u{Vr) 


\/ v{Vr) + 1 


n 


l|A° 


A||l 2 (q^), Vn > 1 


where A minimizes J**. v(Vr) = \\'Pr\\c(L'^(QT))\YPr ^IIc(l2(Qt)) denotes the condition number 
of the operator Vr- As discussed in m Section 4.4], the condition number of Vr restricted to 
A/j C LP‘{Qt) (which coincides with the condition number of the matrix B^Ar^hB]^) behaves 
asymptotically as CT^h~^, where Cr is the constant appearing in (54). This quadratic behavior 
is the typical one for well-posed elliptic problems. Table reports the number of iterations of 
the algorithm, initialized with A° = 0 in Qt- We take e = 10“^° as a stopping threshold for 
the algorithm (the algorithm is stopped as soon as the norm of the residue (/" given here by Ly" 


satishes ||5 "||l2(Qt) — Id ^°||( 7 °||l 2 (q,j,)). We observe that the number of iterates is sub-linear 


with respect to h, precisely, with respect to the dimension nrih = card{{\h\) of the approximated 
problems. This renders this method very attractive from a numerical point of view. 

Tablej^reports the results for r = 1. We get a slightly better estimate for the norm \\Lyh\\L'^(Qr.) 
but slightly worst estimate for the norm \\y — yh\\L'^{QT)- On the other hand, since r acts as an 
augmentation parameter, the convergence of the conjugate gradient algorithm is faster. We also 
check in view of Remark (specifically, (45) and (54)) that the convergence of Xh toward zero is 
much lower for r = 1 than for r = h^: 


Wy-yhh^QT) _ i.02n \\d^{y - yh)\\L^rT) _ rq(,o.i6s 

Mlhq.) ^ ’’ \\d.y\\LHT.) ^ 

Uhh^iQr) = WLynhHQr) = 


h 

7.07 X 10-2 

3.53 X 10-2 

1.76 X 10-2 

8.83 X 10-3 

4.42 X 10-3 

\\y~yh\\L,2i^Q^-^ 

\\d^(y-yh)\\r'2TT) 

1.63 X 10-2 

7.67 X 10-3 

0.937 

7.74 X 10-3 

861 

57 

6.63 X 10-3 

4.95 X 10-3 

1.204 

3.74 X 10-3 

3 321 

103 

2.78 X 10-3 

3.24 X 10-3 

1.496 

1.72 X 10-3 
13 041 

172 

1.29 X 10-3 

2.16 X 10-3 

1.798 

7.90 X 10-'‘ 
51 681 

337 

5.72 X 10-^ 

1.48 X 10-3 

2.135 

3.60 X 10-4 

205 761 

591 

ll^-!/lli,2(rT) 

\\Byh\\L^(Q.r) 

I|A/i||l2(q,j,) 

card({A?i}) 

H CG iterates 


Table 2: Example EXl - BFS element - r = h? - T = 2. 


h 

7.07 X 10-2 

3.53 X 10-2 

1.76 X 10-2 

8.83 X 10-3 

4.42 X 10-3 

\\y~yh\\L,2i^Q^-^ 

2.25 X 10-2 

1.07 X 10-2 

5.23 X 10-3 

2.62 X 10-3 

1.29 X 10-3 

\\d^(y-yh)\\r'2TT) 

ll^-!/lli,2(rT) 

3.7 X 10-2 

3.46 X 10-2 

3.14 X 10-2 

2.76 X 10-2 

2.37 X 10-2 

\\Byh\\L^(QT) 

0.242 

0.207 

0.171 

0.139 

0.112 

A/i l2(q,j,) 

0.147 

0.142 

0.131 

0.116 

0.101 

H CG iterates 

35 

60 

106 

179 

312 


Table 3: Example EXl - BFS element - r = 1. - T = 2. 

We now discuss the results obtained with the reduced HCT finite element on regular but non- 
uniform triangulations of the rectangle XI x (0,T). Precisely, we consider 5 levels of meshes of Qt 
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t X 

Figure 2: The observation y,^^obs on {1} x (0,T) associated to initial data EXl. Example of mesh 
of the domain Qt- 


Mesh number 

1 

2 

3 

4 

5 

tt elements 
tt points 
h 

688 

375 

7.62 X 10-2 

2 752 

1 437 

3.81 X 10-2 

11 008 

5 625 

1.91 X 10-2 

44 032 

22 257 
9.53 X 10-3 

176 128 

88 545 

4.77 X 10-3 


Table 4: Example EXl - HCT element - Information concerning the meshes of the domain Qt- 


described in Table For each of these meshes, we compute h as the maximum of the diameters 
of the triangles composing the triangulation. The coarsest of this meshes is displayed in Figure 
Tablej^collects the numerical results on the reconstruction of the solution y from the observation 
yv,obs obtained with again r = . We observe a slightly super linear convergence for the variable 

yh and A^,: 



Wy-yhh^QT 

\\y\\L^QT) 

- = C>(/l3-22)^ 


= 


h 

7.62 X 10-2 

3.81 X 10-2 

1.91 X 10-2 

9.53 X 10-3 

4.77 X 10-3 

\\y~Vh\\L2(^Q^'^ 

WyWh'^iQ-p) 

3.67 X 10-2 

1.35 X 10-2 

5.99 X 10-3 

2.63 X 10-3 

1.22 X 10-3 

\\^h\\L'^{QT) 

2.12 X 10-2 

1.08 X 10-2 

5.45 X 10-3 

2.53 X 10-3 

1.18 X 10-3 

\\Lyh\\L^(QT) 

1.88 

2.51 

3.26 

4.13 

5.13 

Kh 

2.15 X 10® 

1.11 X 10^ 

1.03 X 10® 

8.67 X 10® 

6.94 X 10® 


Table 5: Example EXl - HCT element - r = h? - T = 2. 


Figure 1^ depicts the exact solution y computed by 0 and its approximation yu (computed 
with the mesh jj 3). Figure 0 represents the relative error \\y — yh\\L'^{QT)/\\y\\L'^(QT) reported also 
in Table[^ Tablesandalso report the value of the condition number Kh of the matrix Arsh¬ 
in both situations, as it is usual for elliptic problems, Kh behaves quadratically with respect to 
h-\ 

We also emphasize that this variational method which requires a finite element discretization of 
the time-space Qt is particularly well-adapted to mesh optimization. Still for the example EXl, 
Figure depicts a sequence of four meshes of Qt = (0,1) x (0,T): the sequence is initialized 
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Figure 3: Example EXl. Exact solution y and approximated solution yh on the mesh ft 3. 
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Figure 4: Example EXl -T = 2 - Relative error \\y — yh\\L^(Q.r)/\\y\\L^{QT) w.r.t. h for the BFS 
element with r = (+) and r = 1 (*), the HCT element with r = h? (□) and r = 1 (o) 

with the coarsest mesh described in Table which is locally refined near boundary Ft (where the 
observation yu,obs is localized) and near VL x {0} (for a better representation of the initial data). 
The three other meshes are successively obtained by local refinement based on the norm of the 
gradient of yh on each triangle of 7/t. As expected, the refinement is concentrated around the 
lines of singularity of yh traveling in Qt, generated by the singularity of the initial position y^. 
Some information concerning these meshes and the approximation errors obtained where this mesh 
adaptation strategy is employed are reported in Table 


Mesh number 

1 

2 

3 

4 

tl elements 

3 452 

6 262 

16 440 

49 499 

jl points 

1 986 

3 442 

8 629 

25 318 

l|y-yhlli,2{Q,j,) 

1.35 X 10-2 

8.02 X 10-® 

5.5 X 10-® 

4.51 X 10-® 

\\y\\L^(Qj.) 

\\^h\\L^(QT) 

2.11 X 10-2 

1.38 X 10-2 

8.57 X 10-® 

4.56 X 10-® 

\\Lyh\\mQT) 

7.03 

6.98 

7.82 

9.01 

Kh 

7.9 X 10® 

3.05 X 10^ 

4.38 X 10i° 

6.19 X lO^i 


Table 6: Example EXl- HCT element - r = h^-T = 2- Iterative refinement. 
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Figure 5: Example EXl - reduced HCT finite element: Iterative refinement of the triangular mesh 
over Qt with respect to the variable yh 


We end this section with some numerical results for the stabilized mixed formulation ( |55[ ). The 
main difference is that the multiplier A is now approximated in a richer space (see @) leading to 
larger linear systems. Table considers the case of the example EXl with T = 2 and a = 1/2. 
In order to compare with the formulation (36), we take again r = h^. We observe the convergence 
w.r.t. h and obtain similar rates and constants to the ones in Table in particular, we have 

\\y-yh\\L^(Q^) 


\\y\\L'^{Q-i 




Finally, we also check - in contrast with the mixed formulation - that the positive parameter r 
does not affect the numerical results. 


5.2 Reconstruction of the solution - Two dimensional case {N = 2) 

In this section we illustrate the method introduced in Section on a two-dimensional example. 
The procedure is similar but a bit more involved on a computational point of view, since Qt is 
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h 

7.62 X 10-2 

3.81 X 10-2 

1.91 X 10-2 

9.53 X 10-3 

4.77 X 10-3 

\\y~Vh\\L2(^Q^'^ 

WyWh’^iQq,) 

\\A\\l^{Qt) 

\\Lyh\\L^QT) 

Kh 

4.1 X 10-2 

6.96 X 10-3 
1.47 

5.05 X 108 

1.55 X 10-2 

3.67 X 10-3 

1.86 

3.01 X 10® 

6.88 X 10-3 

1.71 X 10-3 
2.28 

2.59 X 10^ 

3.03 X 10-3 

8.42 X 10-^ 
2.76 

1.82 X lO^o 

1.39 X 10-3 

3.98 X 10-'‘ 
3.23 

1.31 X 10“ 


Table 7: Example EXl - HCT element - r = h‘^- T = 2- 


Stabilized mixed formulation (55). 


now a subset of We take again c := 1 in fl and d := 0 in Qt- 

In order to approximate the mixed-formulation (§, we consider a mesh 7^ of the domain 
Qx = n X (0, T) formed by triangular prisms. This mesh is obtained by extrapolating along the 
time axis a triangulation of the spatial domain LI. For all the simulations considered in this section, 
LI is the so called Bunimovich’s stadium (see [7]) and T = 3. Figure |^Left displays the domain LI 
and the part T of the boundary on which the observation is available for t G (0, T) while Figure 
[fright displays an example of mesh of domain Qx- 



Figure 6: Bunimovich’s stadium and the subset T of dLl on which the observation is available 
(Left). Example of mesh of the domain Qx (Right). 




Let Zh be the hnite dimensional space defined as follows 




Th = ')pixi,X2)0{t) G C^iQx) 
iph = 0 on Ex 


for every K = x Kt gTh- 


(72) 


V{KxiX 2 ) is the space of functions corresponding to the reduced Hsieh-Clough-Tocher (HCT) C^- 
element recalled in Section 4.1.1; Q(Kt) is a space of degree three polynomials on the interval Kt 
of the form [tjUj+i] defined uniquely by their value and the value of their first derivative at the 
point tj and tj+i. In other words, Yh is the finite element space obtained as a tensorial product 
between the reduced HCT finite element and the cubic Hermite finite element. We check that on 
each element K = Kx^x 2 ^ Kt, the function Lph is determined uniquely in term of the values of 
Ek ■■= {‘p{ai),ipxi{ai),ipx2{ai),(pt{ai),ipxi,t{0‘i)^Tx2A°'i)u = I,-” ,6} at the six nodes at of K. 
Therefore, dimS^ = 36. 

We consider meshes formed by triangular prisms of the domain Qx = H x (0, T). An example of 
such a mesh associated to the domain Qx is displayed in Figure [fright. This mesh is composed by 
10 261 nodes distributed in 18 060 prismatic elements (this mesh corresponds to the mesh number 
2 described in Table [^. 
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We consider three levels of meshes of the domain Qt formed by the number of prisms and 
containing the number of nodes reported in Table 


Mesh number 

1 

2 

3 

Number of elements 

1 860 

18 060 

158 280 

Number of nodes 

1 216 

10 261 

84 241 

Ax 

1.82 X 10-1 

8.2 X 10-2 

3.95 X 10-2 

At (Height of elements) 

0.2 

0.1 

0.05 

h 

2.7 X 10-1 

1.29 X 10-1 

6.37 X 10-2 


Table 8: Characteristics of the three meshes associated with Qt- 


Comparing to the one dimensional situation described in Section |5.1[ the eigenfunctions and 
eigenvectors of the Dirichlet Laplace operator defined on ft are not explicitly available. Conse¬ 
quently, from a given set of initial data, we define as the ’’exact” solution and note y the solution 
obtained numerically with a very fine discretization, from which we can extract an observation on 
Tt- Precisely, we solve the hyperbolic equation Q using a standard time-marching method: we 
employ a HCT finite elements method in space coupled with a Newmark unconditionally stable 
scheme for the time discretization. 

Here, we solve the hyperbolic equation on the spatial mesh which was extrapolated in time in 
order to obtain the mesh number 3 of Qt- This two-dimensional mesh contains 1 381 nodes and 
2 638 triangles and corresponds to the value Acc « 3.95 x 10“^. As for the time discretization, 
we use the value At = 10“^. We denote y the solution obtained in this way for the initial data 
(yo,yi) e X L‘^{n) given by 


(EX2) 


f -Ayo = 10, in fi 
1 yo = 0, on dfi, 


1 / 1=0 in fi. 


(73) 


From ?/, we then generate the observation y^^obs as the restriction of d^y on Ft- The geometry 
used here allows to compute easily the normal derivative of y on the boundary. Finally, from this 
observation we reconstruct yh as the solution of the mixed formulation (28) on each of the three 
meshes described in Table For this simulations we take the augmentation parameter r = . 

Table displays some norms of yh and \h obtained for the three meshes and illustrates again the 
convergence of the method. 


Mesh number 

1 

2 

3 

\\y~yh\\ 

3.75 X 10-2 

1.53 X 10-2 

1.39 X 10-2 

nyWh'^iQ-j,) 

\\Lyh\\L-^{QT) 

2.07 

1.39 

1.09 


9.3 X 10-® 

4.89 X 10-® 

4.57 X 10-® 


5.9 X 10-2 

2.23 X 10-2 

1.08 X 10-2 


Table 9: Initial data (?/o,2/i) given by (73). 


Figure [^Left displays the initial position yo of (73) while Figure [^Right displays the initial 
position yh(-,0) corresponding to restriction at time t = 0 of the solution yh of the inverse problem. 
Figures correspond to the mesh number 2.. The errors between these two functions are given in 
the last row of Table IH 
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5.3 Reconstruction of the solution and the source /i C H ^(f2) - One 
dimensional case N = 1 


We now consider the reconstruction of both the state and the source from a partial observation, 
as discussed in Section In order to construct explicit solution, we hrst recall that the solution 
of 0 with zero initial condition can be expanded as follows : 

y{x, i) = ^ bp{t) sin(p7ra;) 

p>0 

bpit) := — [ sin(pTT{t — s))fp{s)ds, fp(s) := 2a{s) [ sm(pTrx)iJ,(x)dx. 

Ptt Jo Jn 

In the following examples, we take a{t) = 1 +1 defined on [0, T], T = 2 and Ft = {1} x (0, T). 

We now report the resolution of the discrete mixed formulation ( |6l| ). We first consider a rather 
smooth case, with p G Il^(fl) given by 


(EX3) — - l[o,6(](a;) + ^ l[6(,i](a;), d — 1/3. 

Table [^reports the main norms with respect to h. Concerning the augmentation parameter r, we 
use r = h'^ which leads to slightly better approximation of the function p than r = hJ. We check 
the convergence of the approximations (%jh,Ph) as h tends to 0. In particular, we get 


Wy-ytih^QT) 


IIm ~ 




For instance, for h = 8.83 x 10“^ (precisely, Nx = At = 1/160) leading to a linear system 
with 258 566 unknowns, we get a relative error for p^ of the order of 10“^. This allows a very 
good reconstruction of the corresponding solution y. Figure [^Left depicts the function p and 
its corresponding approximation ph- Figure Right depicts along Lt the function (—A)“^(^ — 
dh)/\\{—A)~^p\\ of magnitude 10“®. Here, A denotes the Dirichlet Laplacian. 

Remarkably, the approach also provides good reconstruction of the solution y when the ob¬ 
servation yp^obs is obtained from less regular p function. We consider the following rather stiff 
examples, respectively in LJ{Vl) and H~^{Ll)-. 


(EX4) p{x) = l[a,b]{x) G L^{{),1). a = 0.2, 6 = 0.5 
(EX5) p{x) = l/y/^GH-\t)A)- 


The corresponding normal derivatives yx{I,-) in II^{0,T; L‘^{dlt)) (see [57]) and in L^(0,T) re¬ 
spectively, are depicted on Figure]^ 
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h 

7.07 X 10-2 3.53 X IQ-^ 1.72 x lO-^ 8.83 x lO-^ 7.07 x lO-^ 

\\y~Vh\\L2(^Q^'^ 

WvWh'^iQ'jp) 

\\du{y-yh)\\L2^r^) 

1.72 X 10-3 5.06 X 10-^ 1.28 x 10-* 3.45 x 10-^ 2.14 x 10-^ 

5.9 X 10-3 1.63 X 10-3 8.3 x IQ-^ 3.79 x 10-^ 1.68 x 10-* 

6.5 X 10-* 2.03 X 10-4 5.17 x 10-^ 1.38 x 10-^ 8.89 x 10-® 

7.44 X 10-2 5.26 x 10-2 3.68 x 10-2 2.63 x 10-2 2.35 x 10-2 


Table 10: Example EX3 - BPS element - r = - T = 2. 




Figure 8: EX3 -Left: Function fjL (full line) and its approximation /r/j (dotted line) along 11; Right: 


(-A) ifi-tih) 

||(-A)-VII„i 




along n. 


Table 11 reports some norms with respect to h for the example (EX4). 
convergence as h —>■ 0 with as expected a lower rate : precisely, we compute 


We observe the 


\\y yh\\L^{QT) _ IIm yh\\H-^(n) _ Q^f^o.95^ 

\\y\\L^(QT) ’ 

For h = 8.83 x 10“^, Figure [TT] depicts the functions n and leading to a relative error equal 
to 2.16 X 10“^. The approximation oscillates around p, and suggests that the convergence (in 
agreement with Theorem 3.1) may not hold pointwise but in a weaker (average) sense. Again, this 
weak convergence of the source term is enough to reconstruct with robustness the solution y. The 
function —A“^(/r — yh)/\\ — of magnitude 10“^ suggests the efficiency of the method 

to reconstruct the spatial term /r of the source /. 


h 

7.07 X 10-2 

3.53 X 10-2 

1.72 X 10-2 

8.83 X 10-3 

7.07 X 10-3 


4.72 X 10-3 

2.34 X 10-3 

5.58 X 10-4 

2.96 X 10-4 

2.18 X 10-4 

Wy^L'^iQ^,) 

IlM-Mfellff-qr!) 

1.53 X 10-* 

7.88 X 10-2 

2.5 X 10-2 

2.16 X 10-2 

1.76 X 10-2 

I|a,.(!/-yh)|li2(r,j,) 

2.76 X 10-3 

1.08 X 10-3 

3.57 X 10-4 

1.4 X 10-4 

1.01 X 10-4 

P'/l Z,2 (Qt) 

6.81 X 10-2 

5.07 X 10-2 

3.4 X 10-2 

2.56 X 10-2 

2.31 X 10-2 


Table 11: Example EX4 - BFS element - r = h'^ - T = 2. 

Similar remarks can be made for the example EX5: we refer to Table [T^ and Figure [TT| 
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Figure 9: Normal derivative d^y\q^ = on (0,T) associated to y{x) = l(o. 2 ,o. 5 )( 2 ;) (left) and 

fj.{x) = Ijy/x (right). 




Figure 10: EX4 -Left: 


Function fi (blue) 




in n. 


and its approximation (red) along Right: 


6 Concluding remarks and perspectives 

The mixed formulations we have introduced here in order to address inverse problems for hyperbolic 
equations seems original. These formulations are nothing else than the Euler systems associated 
to least-squares type functionals and depend on both the state to be reconstruct and a Lagrange 
multiplier. This Lagrange multiplier is introduced to take into account the state constraint Ly—f = 
0 and turns out to be the controlled solution of a hyperbolic equation with the source term {d,yy — 
yv,obs') Ipt’" This approach, recently used m a controllability context in m, leads to a variational 
problem defined over time-space functional Hilbert spaces, without distinction between the time 
and the space variable. The main ingredients are, first a unique continuation type property for the 
hyperbolic equation (assuming some geometric conditions on the measurement zone) allowing to 
prove the well-posedness of the mixed formulation, and second, a (strong) generalized observability 
inequality, allowing to quantify the global reconstruction of the solution. 

At the practical level, the discrete mixed time-space formulation is solved in a systematic way 
in the framework of the finite element theory. The approximation is conformal allowing to obtain 
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h 

7.07 X 10-2 3.53 X 1.72 x lO-^ 8.83 x lO-^ 7.07 x lO-^ 

\\y~Vh\\L2(^Q^'^ 

WvWh'^iQ'jp) 

\\du{y-yh)\\L2^r^) 

1.82x10-2 7.74 x 10-2 3.18 x 10-3 1.87x10-3 1.17x10-3 

31.44 11.27 3.96 1.42 1.02 

2.49 X 10-1 2.88 x IQ-i 1.01 x lO-i 5.68 x 10-2 4.72 x 10-2 

2.99 X 10-1 2.35 x IQ-i 1.91 x lO-i 1.62 x lO-i 1.52 x lO-i 


Table 12: Example EX5 - BPS element - r = - T = 2. 



Figure 11: EX5 -Left: 


(-A)- 


ffi(n) 


Function /r (blue) and its approximation (red) along O; Right: 


the strong convergence of the approximation as the discretization parameters tends to zero. In 
particular, we emphasize that there is no need, contrarily to the classical approach, to prove some 
uniform discrete observability inequality: we simply use the observability inequality on the finite 
dimensional discrete space. The resolution amounts to solve a sparse symmetric linear system : 
the corresponding matrix can be preconditioned if necessary, and may be computed once for all 
as it does not depend on the observation yv,obs- Eventually, the space-time discretization of the 
domain allows an adaptation of the mesh so as to reduce the computational cost and capture the 
main features of the solutions. Similarly, this space-time formulation is very appropriate to the 
non-cylindrical situation. 

In agreement with the theoretical convergence, the numerical experiments reported here display 
a very good behavior and robustness of the approach: the reconstructed approximate solution 
converges strongly to the solution of the hyperbolic equation associated to the available observation. 
Remark that from the continuous dependence of the solution with respect to the observation, the 
method is robust with respect to the possible noise on the data. 

Eventually, since the mixed formulations rely essentially on a generalized observability inequal¬ 
ity, it may be employed to any other observable systems for which such property is available : we 
mention notably the parabolic case usually badly conditioned - in view of regularization property - 
and for which direct and robust methods are certainly very advantageous. This is in contrast with 
the Luenberger approach (mentioned in the introduction) which assume the reversibility in time 
of the equation We also mention that this kind of approach may be used to reconstruct potential 
and coefficient. 
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